MS: Spatiotemporal variation in mountain stream metabolism and nitrogen cycling across contrasting flow regimes

Pause to make a column for biomass to chla ratios And also make biomass in in micro grams to more comparable to chl-a

Additive columns for combined surface and pore water nutrient as tot_ in ugL

Kernal density plots of GPP and ER

Colored by years

Modified figure 2

Focal Question:

How does stream metabolism vary based on antecedent flow conditions, and nitrogen availability across small and large stream?

Can we frame wet years as a shrinking productivity window? Or is it more a function of a mismatch in phenology ques, where peak light occurs at unsuitable streamflows (e.g. to high) for epilithic algae?

Ideas:
  1. Figure out what flow threshold at each site we can’t really model GPP at?

    1. How large is that window from year to year?

    2. How is the cumulative GPP during that time?

    3. How is the cumulative PAR and water temp?

  1. Center dynamics on that flow threshold
    1. X-axis could be time since peak “snow pack metric” or “flow metric”
  1. Double check there isn’t evidence of UV suppression, should double check magnitude of surface PAR

Data resolution for DO and metab

DO in black, metabolism in green

### big grid 
WY_n_grid2 <- ggarrange(
  templt,
  templt_BWU,
  templt_GBL,
  templt_GBU,
  ncol = 1, nrow = 4)

WY_n_grid2

Get snotel data

Rough time series plots to start getting at the antecedent flow conditions

(1) Figure out what flow threshold at each site we can’t really model GPP at?

flow stats for when GPP was modeled

Max Q where GPP was obtained…

  • BWL : 2.763 cms

  • BWU : 5.990 cms

  • GBL : 0.622 cms

  • GBU : 0.299 cms

Get DOY dates for max light and flow

Manually identify baseflow onset

Looking at Cumulative GPP

Dotted lines for onset of baseflow Dashed lines for peak PAR

Looking at cumulative GPP + light + flow at each reach:

BW Lower

BW Upper

GB Lower

GB Upper

Trying to parse out why the GPP signal isn’t always decernamble:

Logistic flow plots

Where 0 = DO/Temp data but either no metabolism estimate or GPP/ER = 0

Where 1 = DO/Temp data and GPP > 0 and ER < 0

GPP and flow thresholds?

BWL_logQ_plt <- ggplot(BWL_logQ_df_gpp %>% filter(water_year < 2025), aes(x = log(Q_m), y = GPP_bi)) +
  # Points colored by GPP_mean
  geom_point(aes(fill = jday), color = "black", shape = 21, alpha = 0.9) +  
  #scale_fill_gradient(low = "#b3d4a9", high = "#284a31") +  # Gradient for GPP_mean
  scale_fill_viridis_c(name="day of year") +  # Gradient for GPP_mean
  
  # Logistic regression curve colored by water_year
  geom_smooth(aes(color = as.factor(water_year)), method = "glm", method.args = list(family = "binomial"), se = FALSE, linewidth=0.5) +
  scale_color_manual(values = c("2021" = "#D62828", "2022" = "#F77F00", "2023" = "#003049", "2024" = "#d1ac00")) + # Custom colors for years

  # Labels and theme
  ylab(expression(GPP[bi]~(0~or~1))) +
  xlab(expression(log(Q[m]))) +
  theme_classic() +
  #facet_grid(water_year ~ ., scales = "free") +
  guides(fill = guide_colorbar(title = "Day of Year"), 
         color = guide_legend(title = "Water Year"))  +
  labs(caption = "BW Lower")

Logistic growth lines color by year, individial points color by day of year (light= later)/

ER and flow thresholds?

Table of flow thresholds where metabolism is unreasonable to model:
# Function to fit logistic regression and extract parameters
fit_logistic_stats <- function(data) {
  model <- glm(GPP_bi ~ log(Q_m), data = data, family = binomial)
  coef_vals <- coef(model)  # Extract coefficients
  
  # Calculate steepness (slope) and midpoint
  slope <- coef_vals[2]  
  midpoint <- -coef_vals[1] / coef_vals[2]  

  return(data.frame(
    water_year = unique(data$water_year),
    slope = slope,
    midpoint = midpoint
  ))
}

BW Lower

GPP:

##              water_year      slope   midpoint
## log(Q_m)...1       2021 -2.0554356 -0.2385847
## log(Q_m)...2       2022 -1.1564661 -2.3385433
## log(Q_m)...3       2023 -0.2715363 -4.9852004
## log(Q_m)...4       2024 -1.2827351 -1.8693248

ER:

##              water_year      slope   midpoint
## log(Q_m)...1       2021 -2.0554356 -0.2385847
## log(Q_m)...2       2022 -1.1564661 -2.3385433
## log(Q_m)...3       2023 -0.2715363 -4.9852004
## log(Q_m)...4       2024 -1.2827351 -1.8693248

BW Upper

GPP:

##              water_year       slope    midpoint
## log(Q_m)...1       2021   3.5374237 -4.26012954
## log(Q_m)...2       2022  -1.0596443 -4.11985860
## log(Q_m)...3       2023   0.5265648 -0.06121537
## log(Q_m)...4       2024 -12.0611404 -3.29535264

ER:

##              water_year       slope   midpoint
## log(Q_m)...1       2021   3.5374237 -4.2601295
## log(Q_m)...2       2022  -1.0596443 -4.1198586
## log(Q_m)...3       2023   0.3988148 -0.9032967
## log(Q_m)...4       2024 -12.0611404 -3.2953526

GB Lower

GPP:

##              water_year       slope   midpoint
## log(Q_m)...1       2021  0.34117918   0.516245
## log(Q_m)...2       2022  0.49547049  -1.821918
## log(Q_m)...3       2023 -0.03485753 -36.675012
## log(Q_m)...4       2024 -0.15859985  -6.455164

ER:

##              water_year       slope   midpoint
## log(Q_m)...1       2021  0.34117918   0.516245
## log(Q_m)...2       2022  0.49547049  -1.821918
## log(Q_m)...3       2023 -0.03485753 -36.675012
## log(Q_m)...4       2024 -0.15859985  -6.455164

GB Upper

GPP:

##              water_year      slope  midpoint
## log(Q_m)...1       2021 -2.1348308 -3.430910
## log(Q_m)...2       2022 -0.1903413 -4.478703
## log(Q_m)...3       2023 -0.8599841 -4.445702
## log(Q_m)...4       2024 -0.6843997 -2.921234

ER:

##              water_year      slope  midpoint
## log(Q_m)...1       2021 -2.1348308 -3.430910
## log(Q_m)...2       2022 -0.1898336 -4.421158
## log(Q_m)...3       2023 -0.8599841 -4.445702
## log(Q_m)...4       2024 -0.6843997 -2.921234

Water temp thersholds?

Water temp thersholds?

Flow thresholds based on midpoint in metabolism logistic plots

BWL_logQ_df <- BWL_logQ_df_gpp %>%
  mutate(Q_bi = case_when(
    (water_year == 2021 & Q_m > 1.27) ~ 0,
    (water_year == 2022 & Q_m > 0.10) ~ 0,
    (water_year == 2023 & Q_m > 0.007) ~ 0,
    (water_year == 2024 & Q_m > 0.154) ~ 0,
    (water_year == 2021 & Q_m < 1.27) ~ 1,
    (water_year == 2022 & Q_m < 0.10) ~ 1,
    (water_year == 2023 & Q_m < 0.007) ~ 1,
    (water_year == 2024 & Q_m < 0.154) ~ 1,
    TRUE ~ 0
  ))




BWU_logQ_df <- BWU_logQ_df_gpp %>%
  mutate(Q_bi = case_when(
    (water_year == 2021 & Q_m > 0.015) ~ 0,
    (water_year == 2022 & Q_m > 0.016) ~ 0,
    (water_year == 2023 & Q_m > 0.940) ~ 0,
    (water_year == 2024 & Q_m > 0.037) ~ 0,
    (water_year == 2021 & Q_m < 0.015) ~ 1,
    (water_year == 2022 & Q_m < 0.016) ~ 1,
    (water_year == 2023 & Q_m < 0.940) ~ 1,
    (water_year == 2024 & Q_m < 0.037) ~ 1,
    TRUE ~ 0
  ))



GBL_logQ_df <- GBL_logQ_df_gpp %>%
  mutate(Q_bi = case_when(
    (water_year == 2021 & Q_m > 0.161) ~ 0,
    (water_year == 2022 & Q_m > 0.161) ~ 0,
    (water_year == 2023 & Q_m > 0.161) ~ 0,
    (water_year == 2024 & Q_m > 0.0015) ~ 0,
    (water_year == 2021 & Q_m < 0.161) ~ 1,
    (water_year == 2022 & Q_m < 0.161) ~ 1,
    (water_year == 2023 & Q_m < 0.161) ~ 1,
    (water_year == 2024 & Q_m < 0.0015) ~ 1,
    TRUE ~ 0
  ))


GBU_logQ_df <- GBU_logQ_df_gpp %>%
  mutate(Q_bi = case_when(
    (water_year == 2021 & Q_m >  0.0323) ~ 0,
    (water_year == 2022 & Q_m >  0.0113) ~ 0,
    (water_year == 2023 & Q_m > 0.0117) ~ 0,
    (water_year == 2024 & Q_m > 0.0539) ~ 0,
    (water_year == 2021 & Q_m <  0.0323) ~ 1,
    (water_year == 2022 & Q_m <  0.0113) ~ 1,
    (water_year == 2023 & Q_m < 0.0117) ~ 1,
    (water_year == 2024 & Q_m < 0.0539) ~ 1,
    TRUE ~ 0
  ))

Q ~ K600 relationships

Code sandbox: